local function _gen_erfcinv(erfc)
    local function erfcinv(y)
        if y >= 2 then
            return -100
        elseif y <= 0 then
            return 100
        end
        local zero_point = y < 1
        if not zero_point then
            y = 2 - y
        end
        local t = math.sqrt(-2 * math.log(y / 2.))
        local x = -0.70711 * ((2.30753 + t * 0.27061) / (1 + t * (0.99229 + t * 0.04481)) - t)
        for i = 1, 2 do
            local err = erfc(x) - y
            x = x + err / (1.12837916709551257 * math.exp(-(x ^ 2)) - x * err)
        end
        return zero_point and x or -x
    end
    return erfcinv
end

local function _gen_ppf(erfc)
    local erfcinv = _gen_erfcinv(erfc)
    local function ppf(x, mu, sigma)
        mu = mu or 0
        sigma = sigma or 1
        return mu - sigma * math.sqrt(2) * erfcinv(2 * x)
    end
    return ppf
end

local function erfc(x)
    local z = math.abs(x)
    local t = 1 / (1 + z / 2)
    local r = t * math.exp(-z * z - 1.26551223 + t * (1.00002368 + t * (
        0.37409196 + t * (0.09678418 + t * (-0.18628806 + t * (
            0.27886807 + t * (-1.13520398 + t * (1.48851587 + t * (
                -0.82215223 + t * 0.17087277
            )))
        )))
    )))
    return x < 0 and (2 - r) or r
end

local function cdf(x, mu, sigma)
    mu = mu or 0
    sigma = sigma or 1
    return 0.5 * erfc(-(x - mu) / (sigma * math.sqrt(2)))
end

local function pdf(x, mu, sigma)
    mu = mu or 0
    sigma = sigma or 1
    return (1 / math.sqrt(2 * math.pi) * math.abs(sigma) *
            math.exp(-(((x - mu) / math.abs(sigma)) ^ 2 / 2)))
end

local ppf = _gen_ppf(erfc)

return {
    cdf=cdf,
    pdf=pdf,
    ppf=ppf,
}